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Abstract 



The effect of changing the lattice action for the gluon field on the recently observed 
first order phase transition near zero quark mass is investigated by replacing the Wil- 
son plaquette action by the DB W2 action. The lattice action for quarks is unchanged: 
it is in both cases the original Wilson action. It turns out that Wilson fermions with 
the DBW2 gauge action have a phase structure where the minimal pion mass and the 
jump of the average plaquette are decreased, when compared to Wilson fermions with 
Wilson plaquette action at similar values of the lattice spacing. Taking the DBW2 
gauge action is advantageous also from the point of view of the computational costs of 
numerical simulations. 
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1 Introduction 



A basic feature of the low energy dynamics in Quantum Chromodynamics (QCD) is the 
spontaneous chiral symmetry breaking implying the existence of light pseudo-Goldstone 
(pseudoscalar) bosons. The associated phase structure near zero quark masses has to 
be reproduced in the continuum limit by the lattice-regularized formulations but it 
is in general modified by lattice artifacts at non-vanishing lattice spacing. In lattice 
theories based on Wilson-type quark actions the possible phase structures have been 
investigated up to 0(a 2 ) in the lattice spacing a by Sharpe and Singleton [2] in the 
framework of low-energy chiral Lagrangians [HI E] and using the effective continuum 
description of cut-off effects |H] . Their results allow two possible "scenarios" : the 
existence of the Aoki phase (7j or, alternatively, a first order phase transition near zero 
quark mass. 

In a recent numerical simulation |Hj the phase structure of lattice QCD with 
Wilson fermions and Wilson gauge action has been investigated with the help of the 
twisted mass Wilson fermion formulation For fixed values of a, smaller than 

a ~ 0.2 fm, evidence for a first order phase transition line, near zero quark mass in 
the plane of untwisted and twisted quark mass, has been found corresponding to the 
"second scenario" of ref. • It is important to remark that this line is finite and ends 
at a particular value of the twisted quark mass fj, c . This implies metastability and 
a non-zero minimum of the absolute value of quark- (and pion-) masses. These are 
lattice artifacts which are expected to vanish in the continuum limit where /i c = and 
the first order phase transition line shrinks to a singular point. (For generalizations of 
the results of *Z for non-zero twisted mass see ^2^1E3) Considering, besides the 
bare quark masses, the bare gauge coupling, too, near the continuum limit the first 
order phase transition spans a surface, as it is schematically shown by figure ^ 

It might be speculated that at the microscopic level the occurence of the first 
order phase transition at a > is accompanied by a massive rearrangement of small 
eigenvalues of the Wilson-Dirac operator. The detailed properties and, in particular, 
the strength of the first order phase transition does probably depend on the number 
and distribution of these eigenvalues. It is known that some type of small eigenvalues, 
especially real ones, are associated with small topological dislocations of the gauge field. 
A high probability of these dislocations and of the correspoding small eigenvalues is 
presumably a cut-off effect which can be diminished by an appropriate choice of the 
lattice action. In fact, it is known |14l 1151 fTHl 117] that the small topological dislocations 
can, indeed, be suppressed by taking renormalization group improved (RGI) gauge 
actions as the Iwasaki-action ^H] or the DBW2 action |19j . 

In the present paper we try to answer the question whether the combination of 
RGI gauge actions with the Wilson fermion action does shrink the first order phase 



2 



transition line near zero quark mass. Here we restrict ourselves to the study of the 
DBW2 gauge action which has been successfully applied also in dynamical domain wall 
fermion simulations |2U| . The goal of the present paper is to qualitatively show how 
a change of the gauge action will modify the phase structure. Hence, we do not aim 
here at a high precision study. 

The Iwasaki action is often used in dynamical quark simulations by the CP-PACS 
and JLQCD Collaborations, in particular, in combination with the Sheikholeslami- 
Wohlert clover-improved Wilson fermion action 6 . Earlier results of the JLQCD 
Collaborations indicate |2j that, indeed, a metastability seen in the average plaquette 
can be suppressed by replacing the Wilson plaquette action by the Iwasaki action. (See 
also HU, an( i f° r a review of earlier results on the phase structure of QCD 23 . An 
early discussion of the phase structure of QCD can also be found in |24|.) 

The plan of this paper is as follows: in the next section the lattice action and some 
parameters of the update algorithm are defined. In section |3] we present the results of 
the numerical simulations. Section |I] is devoted to the investigation of the eigenvalue 
spectrum of the Wilson-Dirac operator near the origin. The last section contains some 
discussion and concluding remarks. 



2 Lattice action and simulation algorithm 

2.1 Lattice action 

We apply for quarks the lattice action of Wilson fermions, which can be written as 



f 1 ±4 



iM } ■ (i) 



Here the ("untwisted") bare quark mass in lattice units is denoted by 

H K = am + 4r = — , (2) 
2k 

r is the Wilson-parameter, set in our simulations to r = 1, am® is another convention 
for the bare quark mass in lattice units and k is the conventional hopping parameter. 
In Q the twisted mass \x is also introduced. U XfM G SU(3) is the gauge link variable 
and we also defined U x -^ = U^._- and 7_ M = — 7 M . 

For the SU(3) Yang-Mills gauge field we apply the DBW2 lattice action ^2] which 
belongs to a one-parameter family of actions obtained by renormalization group con- 
siderations. These actions also include, besides the usual (1 x 1) Wilson loop plaquette 
term, planar rectangular (1 x 2) Wilson loops: 

S 9 = PY,\ C ° E {l-^Ret/^l + d ]T {l-^Rel^ 2 }] , (3) 



with the normalization condition Co = 1 — 8c\. (The notation co,i is conventional. Of 
course, c\ should not be confused with the parameter c\ in the effective potential of 
refs. j2m0 The coefficient c\ in (j3J) takes different values for various choices of RGI 
actions, for instance, 

{—0.331 Iwasaki action, 
(4) 
-1.4088 DBW2 action. 

Clearly, c\ = corresponds to the original Wilson gauge action with the plaquette 
term only. Note that for c\ = —1/12 one obtains the tree-level improved action in the 
Symanzik-improvement scheme |25) . 



2.2 Twist angle 

An important quantity is the twist angle u, the polar angle in the plane of the untwisted 
and twisted mass We present here a method which allows to determine the 

twist angle only on the basis of symmetry of the correlators defined in a given point of 
bare parameter space (see also (Hj). 

Following we introduce the twist angle oj as the chiral rotation angle between 
the renormalized (physical) vector and axialvector currents V£„, A x ^ and the bare 
bilinears of the x-fields V£„, A x ^: 

= Xx\ T al^Xx , A a xfi = Xx\ T al^>Xx ■ (5) 

With the renormalization constants Zy and Z A we have 

V?n = z vV^ cos w + e ab Z A A b Xfl sin u , (6) 

A a x ^ = Z A A a X(i cos u + e ab Z V V^ sin u (7) 

where only charged currents are considered (a=l,2). 

The twist angle to is related to the ratio of the renormalized twisted and untwisted 
masses entering the chiral Ward identities [TO] . (In ^U] this definition of the twist angle 
was called a.) We define, in addition, the two auxiliary angles 

ujy = arctan(Z^^y 1 tanw) , uja = arctan(ZyZ^ 1 tanw) . (8) 

In terms of uy, oja eqs. © and Q are written as 

= My (cos ujy + e ab sin uj v A b XfM ) , (9) 
A% = M A (cos uaA^ + e ab sin u A V^) (10) 

where the overall multiplicative renormalization is (X = V, A) : 



From (JHJ it follows that 

uj = arctan {^J tan ujy tan uja^J ■ (12) 

As shown by the relations in Q and (|12j) . the values of cj, and uja coincide for 
\uj\ = 0, 7r/2. However, for other angles they are, in general, different and the difference 
goes to zero in the continuum limit only as fast as Zy /Za — * 1. 

A possibility to determine ujy and uja is to impose the vector and axialvector Ward 
identities, respectively, with a suitable insertion operator O x . For instance, in the 
vector case one can use the Ward identity 

\d u VlT, 0... ) = => tanwy = - . . (13) 

x,y l^x^^O^xO I 

Here the indices + and — refer to the charged components t± = t\ ± 1T2 and <9* denotes 
the backward lattice derivative. 

Another possibility for determining the twist angles toy, uja and to is to impose 
parity conservation for suitable matrix elements, for instance with the pseudoscalar 
density P± = Xx^lbXx- 

zZ( A toVyo) = E^JoPy) = 0. (14) 

x,y x,y 

These equations admit the solution 

tanUv = -^jiV+oPj) ? (15) 
Y,x,y\A x o Py ) 

tanuJA = *W +tanuJ v^,y(At A~ ) ^ 

In ()14|) one can also take the derivatives of the currents instead of the currents 
themselves. For instance, taking the divergence of the vector current in the second 
equality gives the same equations as ((T3*|) with O = P. 

Once ujy and uja are determined, the twist angle uj can be obtained by eq. (|P2*|) . 
This method for determining the twist angle can also be used in case of simulations 
with partially quenched twisted mass quarks. The estimate of uj is, of course, affected 
by 0(a) ambiguities. For non-zero twisted mass the critical bare untwisted quark 
mass fi K = [i KCT -, or the critical hopping parameter k ct = (2fi Kcr )~ 1 , is signaled by 
\uj\=tt/2. 



2.3 Updating algorithm 

Concerning updating in our numerical simulations, we apply the two-step multi-boson 
(TSMB) algorithm which has been tuned to QCD simulations with Wilson quarks 
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in previous works |23 EHl EH1 EDI E] • (For details and references see in these papers.) 
In ref. [27] there is an approximate formula for the computational cost of an update 
cycle in terms of Matrix- Vector-Multiplications (MVMs): 

N M VM/cycle ~ 6(rt B rtiiV$ + Nu) + 2n B (n 2 + n 3 )iVc + IgFg ■ (17) 

Here n\ 23 are the orders of polynomials used in the two approximation steps, tib gives 
the multiplicity in determinant breakup, N$ is the number of local bosonic sweeps 
per update cycle, Nu the number of local gauge sweeps, Nc the number of global 
Metropolis accept-reject correction steps, and Iq and Fq give the number of MVMs in 
the global boson heatbath and its frequency, respectively. 

The number of MVMs can also be converted into the number of floating point 
operations by noting that in our code, for vanishing twisted mass, we have 

1 MVM ~ 1.2 • 10 3 Offop , (18) 

where Q is the number of lattice points. For non-zero twisted mass there is an additional 
factor 2 due to the flavour index. (This does, however, not mean that twisted mass 
fermions are a factor of two more expensive since in this case the two flavours are 
incorporated in one fermion matrix and the polynomial approximations have lower 
orders: see appendix A. 2 of p.) 

Measuring the integrated autocorrelations Ti n t as a function of the quark mass in 
lattice units am q and of the lattice volume O, previous experience tells that one can 
approximate the computational cost of a number of update cycles equal to Ti nt by the 
simple formula 

C Tint ~F(am q )- z n. (19) 

According to ref. jUJ , in case of combining the Wilson fermion action with the Wilson 
plaquette gauge action, the power of the inverse quark mass is close to z = 2. The 
overall factor F depends on the quantity under investigation. For Wilson quarks with 
Wilson gauge action the previous results can be summarized, for instance, for the 
average plaquette and for the pion mass determined with a randomly chosen source by 

M 

Fplag ^ 7 • 10 6 flop , F mw ~ 10 6 flop . (20) 

Let us note that the approximate formula in (|19|) has been, up to now, verified only for 
some fixed values of the gauge coupling j3. The /3-dependence of F has not yet been 
systematically investigated. 

3 Numerical simulation results 

Our aim is to compare the phase structure of two-flavour {Nf = 2) QCD near zero 
quark mass for Wilson lattice fermion action and DBW2 gauge action with the one 
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observed in ref. p] EI for Wilson fermion action and Wilson (plaquette-) gauge action. 
Since the phase structure obviously depends on the lattice spacing, we have to find the 
values of the bare parameters (j3,fjt K ) in the lattice action which correspond to 

quark mass m q ~ and to the same lattice spacing as in P3IBJ, namely a ~ 0.2 fm. For 
having a fair comparison, the lattice volume has to be kept constant, too, because the 
metastability phenomenon does also depend on it. Therefore, we shall compare the 
results on 12 3 x 24 lattices. 

A possibility for facilitating the parameter tuning is to explore the position of 
the high-temperature phase transition on lattices with time extension JVj = 4 and 
Nf = 6 for small quark masses, which mark a = 0.25 — 0.30 fm and a = 0.17 — 0.20 fm, 
respectively. (This method with N-t = 4 has been applied, for instance, in |27|.) A 
useful first orientation is also provided by the quenched studies. (For a useful collection 
of data on RGI gauge actions see [HJ an d references therein). For specifying the actual 
value of the lattice spacing we determine the Sommer scale parameter in lattice units 
ro/a .32], which we set by definition to be ro = 0.5 fm, independently from the quark 
mass. 

In order to localize the Nt = 4 high-temperature phase transition we fixed the 
gauge coupling at f3 = 0.55 and changed the bare quark mass fj, K (or, equivalently, the 
hopping parameter k = (2fi K )~ 1 ). The results on an 8 3 x 4 lattice for the absolute 
value of the Polyakov-line and average plaquette are given in figure As it is shown 
by the figure, the transition with the DBW2 action is rather smooth, barely visible. 
This has to be contrasted with the strong and sudden increase of both Polyakov-line 
and average plaquette in case of the Wilson plaquette action, which is also shown for 
comparison in figure 121 

A similar analysis on 12 3 x 6 lattices at (3 = 0.67 gives qualitatively similar results 
but there the difference between the DBW2 and the Wilson plaquette action is smaller 
because the transition for the Wilson plaquette action becomes weaker. 

3.1 Results on an 8 3 x 16 lattice at (3 = 0.55 

The runs on an 8 3 x 16 lattice at j3 = 0.55 and \i = were started from the low- 
temperature phase by taking four copies in the time direction of some of the 8 3 x 4 
lattices. The parameters of these runs are specified in the first part of table E Be- 
sides the hopping parameter n also some parameters of the TSMB updating algorithm 
are specified: the orders of the polynomials used ni,2,3 and the interval covering the 
eigenvalues of the squared preconditioned hermitean quark matrix [e, A]. 

In the 8 3 x 16 runs we looked for signals of metastability but we did not find any. The 
results for some interesting quantities are collected in the first part of table EJ the pion 
(i.e. pseudoscalar meson) and p-meson masses and the bare quark mass in lattice units 
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am^ CAC . Some of these quantities are also shown in figure El The scale parameter in 
lattice units ro/a was also determined. We note in passing that at this small value of (3 
and with our partly low statistics the evaluation of tq is rather difficult. Nevertheless, 
in order to estimate quantities also in physical units, we performed a purely statistical 
analysis for ro, being aware of the fact that systematic effects can be large. 

The bare quark mass am^ CAC is defined by the PCAC-relation containing the 
axialvector current A" in © and the pseudoscalar density insertion: 

pcac = (d;A^P y ) 

Since for the moment we do not determine the Z-factors of multiplicative renormaliza- 
tion, the bare quark mass am^ CAC contains an unknown 0(1) Z-factor Z q = Zp/Za- 
In the following analysis we extracted the quark mass with the method detailed in 
ref. [27], see section 3.1.1 there. 

In agreement with the absence of a signal for metastability, the //^-dependence of 
the pion mass and quark mass in figure El is consistent with the absence of a first order 
phase transition at this gauge coupling {(3 = 0.55). A rough estimate for the value of 
the lattice spacing is a ~ 0.30 fm in the positive quark mass phase and a ~ 0.23 fm 
in the negative quark mass phase. The upper panel in figure 01 suggests the existence 
of a short interval (/U K £ [2.62,2.63] or k G [0.190,0.191]) of an Aoki phase near zero 
quark- and pion- masses. This behaviour is qualitatively similar to the one for the 
Wilson plaquette action which also shows the existence of the Aoki phase at strong 
gauge coupling [3*3*] . 

3.2 Results on a 12 3 x 24 lattice at (3 = 0.67 

With a short investigation of the high-temperature phase transition on a 12 3 x 6 lattice 
one can easily localize the gauge coupling and bare quark mass /i K = (2k) _1 where 
the lattice spacing is about a factor 3/2 smaller than at (3 = 0.55. It turned out that 
one can take (3 = 0.67 and k ~ 0.17. Fixing (3 = 0.67 and changing k we performed 
several runs on a 12 3 x 24 lattice. In this way the physical volume of the lattice is 
approximately the same as the one of an 8 3 x 16 lattice at (3 = 0.55. In order to be 
able to compare with the results of ref. PP, besides the runs with fi = 0, at this (3 we 
also considered a non- vanishing twisted mass \i = 0.01. 

First we looked also here at \i = for a signal of metastability in the average 
plaquette and we found it near k = 0.167 — 0.168, as it is shown by the upper panel of 
figure |31 Note that the average plaquette values are substantially higher here than at 
/3 = 5.2 with the Wilson plaquette gauge action in pQ: A p i aq = (|ReTr U p i aq ) ~ 0.59 
instead of A p i aq ~ 0.52. This qualitatively shows that the gauge field with DBW2 is 
smoother. 



We also determined the pion, p-meson and quark masses, with the results given 
in table [21 (For a graphical representation of some of these results see also the upper 
panels of figure [SJ) For the extraction of ro/a we performed only a statistical analysis 
also here, neglecting the systematic effects. Let us give a range of values for the runs 
in table EJ For A\ to we find 2.37 < ro/a < 2.76 in the low and for to 
2.72 < ro/a < 3.17 in the high plaquette phases, respectively. For A^B^ we find 
2.39 < r /a < 2.54 in the low and for C h to £ h we find 2.89 < r /a < 3.07 in the high 
plaquette phase. 

From the values of the scale parameter ro/a we determined the lattice spacing, and 
found a ~ 0.18 — 0.21 fm in the positive and a ~ 0.16 — 0.18 fm in the negative quark 
mass phase, respectively. This is quite close to the values obtained in both phases with 
the Wilson plaquette gauge action at (3 = 5.2 in ref. |Tj. 

Going to the positive twisted mass /i = 0.01, the metastability in the average 
plaquette disappears on our 12 3 x 24 lattice, as it is shown by the lower panel of figurellj 
Having in mind the strong metastability signal in the average plaquette observed on 
a 12 3 x 24 lattice at (3 = 5.2 and n = 0.01 with the Wilson plaquette gauge action in 
, the absence of the metastability here signals a dramatic improvement of the phase 
structure due to the DBW2 gauge action. The presence of metastability at // = and 
the absence of it at fj, = 0.01 indicates the existence of a rather short first order phase 
transition line near the origin in the (/i K , /i)-plane. Of course, for a precise localization 
of the first order phase transition line a detailed study of the infinite volume limit is 
required, which is beyond the scope of this paper. 

An important question is the minimal value of the pion mass ra™ m associated to 
the first order phase transition line. A precise definition of m™ n could be the value of 
the infinite volume pion mass just at the position of the first order phase transition, 
defined by the equal depth of the two free energy minima in the infinite volume limit. 
To obtain this would be rather demanding. Although the volume dependence could 
be studied beyond our volume extension of L ~ 2.4 fm, for instance on a 16 3 x 32 
lattice, the precise comparison of the free energy minima would be quite difficult. An 
approximate determination of m™ in can be obtained by requiring the equality of the 
pion mass in lattice units am^ in the two phases on our 12 3 x 24 lattices. For this a 
linear extrapolation of (am^) 2 from the points on both sides of the phase transition 
can be considered. As shown in the upper left panel of figure 03 our result at \x = 
is (am^) 2 = 0.0881. This implies, with the range of ro values given above, that 
m mvn ^ 251 MeV in the positive quark mass phase and m™ m ~ 374 MeV in the phase 
with negative quark mass. 

The minimal charged pion mass at fj, = 0.01 is ra™" ~ 360 MeV (see lower left 
panel of figure . This originates from the non-zero value of the twisted quark mass 
and not from the presence of a first order phase transition. 
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In the right panels of figure El the bare quark mass in lattice units is shown. The 
dashed lines are linear fits to the points with positive and negative quark mass, respec- 
tively. At zero twisted mass (upper panel) the metastability region near the first order 
phase transition is clearly displayed. At fi = 0.01 (lower panel) the difference between 
the two dashed lines is smaller. This difference may be interpreted as a consequence 
of a cross-over in the continuation of the first order phase transition line. In the figure 
there is also a linear fit to all points shown (full line) which goes reasonably close to 
every point. The two dashed lines also give lower and upper bound estimates for the 
critical hopping parameter: 0.1661 < K cr < 0.1689. 

Another way to estimate the critical hopping parameter (i.e. critical bare un- 
twisted quark mass) is to determine the twist angle and find n cr = (2^ Kcr ) _1 where 
it equals ir/2. Considering, for definiteness, the twist angle ujy defined in section [2~2l 
the fit in figure El gives K cr = 0.16651(2), in good agreement with the previous es- 
timate. (Actually the numbers in figure come from the vector Ward identity (|13j) 
but, within errors, (|15|) gives compatible results.) The Z-parameter appearing in this 
fit for uiy is Z Q y = ZaZs/(ZvZp) (see section l2~2l and [TH]). According to figure El 
we have Z Q y = 0.959(30). Since from an analogous fit to uja one could determine 
Z a = ZvZs/(ZaZp), this also offers a relatively easy way to obtain the Z-parameter 
combinations Za/Zv and Zp/Zs- 

The quantities {r^m^) 2 and am^ CAC can also be plotted against each other (see 
figure and El for fi = and fi = 0.01, respectively). Figure and the data in table [3 
show that at fi = 0, in the metastable region beyond the minimal pion mass, one can 
also reach values close to the physical value ~ 140 MeV. 

3.3 Topological charge 

The RGI gauge actions, and in particular the DBW2 action, are known to slow down the 
transitions between different topological sectors both in quenched |16j and in dynamical 
domain wall simulations [20]. In order to check this we determined the topological 
charge Qt p in some of the runs using a cooling method [HI]. In the following we 
denote the result from the cooling analysis by "topological charge" , being aware of the 
fact that this definition contains some degree of arbitrariness. However, for our aim of 
testing the autocorrelation time this definition is sufficient. 

In the run with label (C h ) (12 3 x 24 lattice, f3 = 0.67, \i = 0.01, k = 0.167) the 
history of the topological charge is shown in the upper panel of figure E3 (The lower 
panel is a histogram of Qtop-) The analyzed configurations are separated by 10 TSMB 
update cycles. In this point, according to table the quark mass is about m q ~ 
0.3m s trange and the pion mass ~ 380 MeV. As the figure shows, the topological 
charge is often changed. Its integrated autocorrelation in this run is r*°? ~ 180, but 
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there is obviously a long tail of the autocorrelation which is not yet properly taken into 
account in a run of this length. In any case, r*°f is substantially longer than those of 
the average plaquette ijf^ q ~ 22) or of the pion mass (t™J ~ 6) in table 

In another run, the one with label (Ci) (12 3 x 24 lattice, (3 = 0.67, fi = 0, k = 
0.167), where the quark mass is about m q ~ 0.18 m s t ra nge an d the pion mass ~ 
295 MeV, the general picture is similar to figure EJ The integrated autocorrelation here 
comes out to be r*°f ~ 70, but this value is even less reliable because the run is shorter. 

In spite of these relatively long autocorrelations, it is clear that in a sufficiently 
long run, say of length 1000 r™^ , which would be needed anyway for a good statistics 
on other quantities, the different topological sectors could be properly sampled by the 
TSMB algorithm. Therefore, at these bare parameter values, there is no problem with 
the suppression of the transitions between different topological sectors. 

3.4 Results about the update algorithm 

In this paper we applied the TSMB update algorithm [2Sj. The estimates of the 
autocorrelations in different runs and the cost estimates obtained using (|17|) are given 
in table|31 Since in our relatively short runs the autocorrelations can only be estimated, 
say, within a factor of two, the numbers in table 01 give only a first orientation. 

Qualitatively speaking, the 12 3 x 24 runs with "low-plaquette" (positive quark mass) 
have lower costs than the corresponding runs with "high plaquette" (negative quark 
mass): at the same absolute value of the bare quark mass the runs in the negative 
quark mass phase have in most cases at least by an order of magnitude higher costs 
than those in the positive quark mass phase. The reason of the higher cost at negative 
quark mass is that the smallest eigenvalues fluctuate more frequently to very small 
values. 

There is also a general tendency that the overall factors F decrease for decreasing 
absolute value of the quark mass. In fact, the data on F show that in the small 
quark mass region an inverse quark mass power z = 1 is a better approximation than 
z = 2, which has been observed in previous simulations with the Wilson plaquette 
gauge action (23 EEl IMl EE IH| At // = the overall factor for the average plaquette 
Fpi aq in the parameterization 



turns out to be F p i aq ~ 2 • 10 7 in the positive quark mass phase and F p i aq ~ (2 — 3) • 10 8 
at negative quark mass. The corresponding numbers at \i = 0.01 are between these 
two values. 

Let us note that at the smallest quark masses a final reweighting correction has to 
be applied because the smallest eigenvalues cannot be always kept in the interval of 
polynomial approximations. Sometimes they fluctuate below the lower limit e. 
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4 Eigenvalue spectra 



Looking at the eigenvalue spectra of the Wilson-Dirac fermion matrix Q at small 
(untwisted) quark masses (see, for instance, in section 4 of ref. |27|) it seems plausible 
that near zero quark mass there has to be a massive rearrangement of eigenvalues. 
This is because in the path integral small eigenvalues are strongly supressed by the 
zero of the fermion determinant. At the sign change of the quark mass the eigenvalues 
have to somehow avoid the zero of the determinant at the origin. It is plausible that 
this eigenvalue rearrangement is related to the phase transition at zero quark mass. 

An interesting question is how the behaviour of eigenvalues in the small quark mass 
region is influenced by a non-zero twisted mass term. 

We investigated the eigenvalue spectra by the Arnoldi method on 8 3 x 16 and 
12 3 x 24 lattices in some of the runs listed in tabled Typically 100-200 eigenvalues 
were determined on 10-30 independent gauge field configurations. The parameters of 
the Arnoldi code were set for searching the eigenvalues with the smallest absolute value. 

The results at \i = on an 8 3 x 16 lattice are shown by figure 1101 In the upper 
panels of the figure, where the quark mass is positive, typical "half-moons" filled with 
eigenvalues can be seen, which correspond to the figures in ref. |22j. At negative quark 
mass - in the lower part of figure EH ~ an almost empty segment in the middle of the 
"half-moon" appears. Comparing the two figures at negative quark mass one can also 
see how this segment is gradually emptied during equilibration. 

It is remarkable that even after equilibration there are some real ("zero- mode") 
eigenvalues on the positive axis. Our Arnoldi code did not find in these configurations 
any negative real eigenvalues. In addition, it is quite surprising that, apart from the 
empty segment in the middle, the half-moon-shaped deformation of the eigenvalue 
region observed at small positive quark masses does not disappear for small negative 
quark masses either. 

The effect of a non-zero twisted mass on the eigenvalue spectrum on an 8 3 x 16 
lattice is illustrated by figure 1111 It can be seen that the strip around the real axis 
— li < Im(A) < +n is free from eigenvalues. Let us remark that also in presence of a 
non-zero twisted mass we studied the spectrum of the operator of equation which 
corresponds to the so called "twisted basis". In the "physical basis" [35], for u = 7r/2, 
the spectrum of the Dirac operator lies in a vertical line parallel to the imaginary axis 
and is shifted from the origin by li (exactly as in the continuum). 

Going to larger (3 (smaller lattice spacing) the visible difference is that the "half- 
moons" are straightened and come closer to the origin: see figure IT2l Otherwise most 
qualitative features are unchanged. There is, however, a marked difference in the 
number of real eigenvalues (for li = 0): in the upper panels of figure ITUl there are lots of 
them, whereas at larger (3, in the upper panel of figured! their number is substantially 
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reduced. 

The effect of changing the gauge action can be seen by comparing figure ^] with 
the eigenvalue spectra in case of the Wilson plaquette gauge action at a similar lattice 
spacing a ~ 0.2 fm in figure ITTT1 The fact that in the case of using the Wilson gauge 
action the pion mass is larger than in the case of the DBW2 action is reflected by a 
movement of the "half- moons" farther away from the origin. 

5 Conclusion 

The main conclusion of this paper is that, indeed, exchanging the Wilson plaquette 
gauge action with the (renormalization group improved) DBW2 action shows substan- 
tial effects on the phase structure: We performed a qualitative study of the phase 
structure of lattice QCD by changing the gauge action and compared the Wilson- 
plaquette and DBW2 actions at a lattice spacing a ~ 0.2 fm in the positive quark mass 
phase and a ~ 0.17 fm in the phase with negative quark mass. This means, at (3 = 5.2 
for the Wilson plaquette action and (3 = 0.67 for the DBW2 action. At this comparable 
situation the metastability signaled by the existence of long living states with different 
average plaquette value and quark masses with opposite sign becomes weaker and the 
minimal pion mass and the jump in the average plaquette between the phases with 
positive and negative quark mass decrease. 

For vanishing twisted mass fx = the metastability occurs in the hopping parameter 
range 0.167 < k < 0.168. Going to the twisted mass value \x = 0.01, which is the same 
as in the numerical simulations of ref. the metastability disappears on our 12 3 x 24 

lattices. It might reappear on larger lattices, but our 12 3 x 24 data indicate that the 
jump in the average plaquette is at least by a factor of ten smaller than the one observed 
in P. 

At a lower [3 value (3 = 0.55, which corresponds to lattice spacings a ~ 0.30 fm and 
a ~ 0.23 fm for positive and negative quark mass, respectively, our simulation data 
are consistent with the existence of the Aoki phase. This is similar to the situation 
for (3 < 4.6 in case of the Wilson plaquette action |33| . The schematic picture of the 
suggested phase diagram in the (f3, k, fi) space, both for DBW2 and Wilson plaquette 
gauge actions, is shown by figure ITU 

The minimal pion mass in a stable phase can be estimated from our simulation data 
at (3 = 0.67 and vanishing twisted mass on a 12 3 x 24 lattice to be m™ m ~ 250 MeV in 
the positive quark mass phase and m™ m ~ 375 MeV in the phase with negative quark 
mass. On larger lattices this value is expected to be 10-20% smaller due to the finite 
volume effects which are non-negligible on the 12 3 x 24 lattice, especially in the negative 
quark mass phase where the lattice extension is only L ~ 2.0 fm. At positive twisted 
mass /i = 0.01 the estimate for the minimal charged pion mass is m™ in ~ 360 MeV, 
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a value entirely due to the non-zero twisted mass and not to the first order phase 
transition. 

Besides the pion- and /9-meson masses, at non-vanishing twisted mass, we also 
determined the twist angle ojy as a function of the bare untwisted quark mass fj, K . 
The /^-dependence of ujy can be well described by the expected arctan-function (see 
figure EJ. From the fit one obtains the value of the critical hopping parameter K cr = 
0.16651(2) and an estimate of a combination of Z-factors. 

In some of the simulation runs we also monitored the history of the topological 
charge (see, for instance, figure EJ). Although the autocorrelation of the topological 
charge is markedly longer than those of the average plaquette or of the pion mass, in 
a good statistics run with, say, thousand times the integrated autocorrelation length 
of the pion mass, the different topological sectors could be properly sampled. 

In order to illustrate the rearrangement of the small eigenvalues near the zero quark 
mass phase transition we investigated in some detail the eigenvalue spectrum of the 
non-hermitean fermion matrix defined in (^Q). To our surprise, the transition from 
positive to negative quark mass is signaled in the eigenvalue spectrum by the opening 
up of an almost empty segment in the "half-moon" occupied by the eigenvalues near 
the origin. The introduction of a non-vanishing twisted mass causes the appearance 
of an empty strip [—(J,, +fj] on both sides of the real axis. The effect of larger j3 is 
to straighten the "half-moon" occupied by the small eigenvalues. At the same time 
the small real eigenvalues at zero twisted mass, which are causing the problem of the 
so-called "exceptional gauge configurations" in partially quenched simulations, occur 
much less frequently. 

A welcome side-effect of introducing the RGI gauge action is the speed-up of the 
TSMB update algorithm. (This presumably also applies to other update algorithms, 
but in this paper we only used TSMB.) This can be qualitatively understood by the 
reduction of the probability for small size "dislocations" in the gauge field and for 
the less frequent occurrence of small real eigenvalues. (This is qualitatively similar 
to the conclusions of ref. j^Sj, obtained in another setup.) The computational cost 
as a function of the quark mass can be better approximated in the small quark mass 
region by an inverse power behaviour of only (am q )~ l than by the behaviour (am q )~ 2 
observed previously with the Wilson plaquette action. 

The results of the present paper indicate that the combination of Nf = 2 Wilson- 
quarks with the DBW2 gauge action leads to a phase structure with a weaker first 
order phase transition than Nf = 2 Wilson-quarks with the plaquette gauge action 
at a comparable value of the lattice spacing. For the moment we have no detailed 
information on the dependence of the phase structure on the parameter c\ in the 
gauge action which multiplies the rectangular Wilson loops. It is possible that the 
optimal choice is different from c\ = —1.4088, for instance, c\ = —0.331 for to the 
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Iwasaki action or c\ = —1/12 for the tree-level- improved Symanzik action. The best 
choice of c\ might also be influenced by the positivity problem of improved actions jHJ 
and /or by the convergence rate of lattice perturbation theory [37] . 

An important open question, which remains to be investigated in the future, is the 
/3-dependence of the phase structure for Wilson-type lattice fermions. It is expected 
that closer to the continuum limit the minimal pion mass and the jump in the average 
plaquette become smaller and finally, in the continuum, the first order phase transition 
line in the plane of untwisted and twisted quark mass shrinks to a first order phase tran- 
sition point. The faster this actually happens the better it is for phenomenologically 
relevant numerical QCD simulations with Wilson-type quarks. 

Another important question is, whether the DBW2 gauge action in combination 
with Wilson twisted mass fermions shows a good scaling behaviour. To this end, a 
simulation at a higher value of (3 than the one used here is necessary. Both questions 
mentioned above are presently investigated by our collaboration. 
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Tables 



Table 1: Bare couplings and parameters of the TSMB algorithm in runs 
with the DBW2 gauge action. The determinant breakup multiplicity is 
ub = 4 in all runs. Small letters label runs on 8 3 x 16 lattices at (5 — 0.55 
whereas capital letters stand for runs on 12 3 x 24 lattices at f3 = 0.67. The 
suffix I and h denote "low" and "high" plaquette phase, respectively. Those 
runs with a caligraphic letter are performed with an additional twisted mass 
term (ji = 0.01). The number of analyzed configurations is given in the 
last column. An asterix on these numbers denotes that a few configurations 
have very low f<C 1 ) reweighting factors. The analyzed gauge configura- 
tions are separated by 10 update cycles, except for run (a) where they are 
separated by 100 update cycles. 
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Table 2: Results of runs specified in tabled for different quantities. 
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Table 3: The cost of an update cycle C cyc ie in thousands of MVMs ac- 
cording to eq. and the estimated integrated autocorrelation lengths in 
update cycles obtained from runs specified by tabled The suffix plaq and 
m n refer to the average plaquette and the pion mass, respectively. The last 
two columns give the factors F calculated from eq. Alty) with z = 2. 
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Figures 




Figure 1: The schematic view of the first order phase transition surface in 
the (j3,K,fi) space close to the continuum, limit. (f3=bare gauge coupling, 
K=hopping parameter, ji=bare twisted quark mass, fi K = (2/t) _1 =bare un- 
twisted quark mass.) The crosses mark the second order boundary line of 
the first order phase transition surface. The strong coupling region near 
P = is not shown in this figure. 
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Figure 2: Upper panels: the signals of the N t — 4 non-zero temperature 
transition on an 8 3 x 4 lattice with the DBW2 gauge action. Lower panels: 
the same with Wilson gauge action. Left panels: absolute value of the 
Polyakov line, right panels: average Wilson loop, both as a function of k. 
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Figure 3: Results of the numerical simulation on an 8 3 x 16 lattice at 
j3 = 0.55: upper panel the square of the pion mass (am^) 2 , lower panel 



the PC AC quark mass am 



PCAC 



In the upper panel the dashed lines are 



extrapolations to zero pion mass: at right it is a linear fit of four points, 
at left a straight line connecting two points with small quark mass. 
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Figure 5: Results of the numerical simulation on a 12 3 x 24 lattice at 
f3 = 0.67 as a function of fi^ = (2k)~ 1 : upper panels fi = 0, lower panels 
H = 0.01. Left panels: (am n ) 2 , right panels: the bare PCAC quark mass 
am^ CAC . The straight lines are fits to the points in the positive and neg- 
ative quark mass phase, respectively. The horizontal line in the upper left 
panel shows the estimated value of the minimal pion mass in lattice units. 
The straight lines in the right panels are explained in the text. 
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Figure 6: T/ie dependence of the twist angle uoy on the bare untwisted 
quark mass fi K = at (3 = 0.67, fi = 0.01 on a 12 3 x 24 lattice. 

The fit determines the critical hopping parameter to be n cr = (2/x Kcr ) _1 = 
0.16651(2). 
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Figure 7: The dependence of ' (r m^) 2 on am^ CAC at (5 — 0.67, ji = on a 
12 3 x 24 lattice. The dashed lines are linear fits of the form A\am^ CAC \. 
At right all three points are included in the fit, at left only two of them with 
the smallest quark masses. The fit parameter is for positive and negative 
quark mass A = 23.6 and A = 43.5, respectively. 
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Figure 8: The dependence of (rom^) 2 on am^ CAC at (3 = 0.67, /i = 0.01 
on a 12 3 x 24 lattice. The dashed line is a spline interpolation for guiding 
the eyes. 
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Figure 9: The topological charge on a 12 3 x24 lattice at (5 — 0.67, ji = 0.01, 
k = 0.167. Upper panel: the run history, with configurations separated by 
10 TSMB update cycles. Lower panel: the distribution of the topological 
charge in this run. 
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Figure 10: Eigenvalues of the Wilson-Dirac fermion matrix (CP with small 
absolute value for fi = 0, (3 = 0.55 on an 8 3 x 16 lattice. Upper left panel: 
k = 0.184. Upper right pannel: k = 0.190. Lower panels: k = 0.194 at the 
beginning of equilibration (left panel) and after equilibration (right panel). 
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Figure 11: Eigenvalues of the Wilson-Dirac fermion matrix (QP with small 
absolute value for fi = 0.01, j3 = 0.55 on an 8 3 x 16 lattice. Upper panel: 
K = 0.184, lower pannel: k = 0.186. 
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Figure 12: Eigenvalues of the Wilson-Dirac fermion matrix with small 
absolute value at (3 = 0.67, k = 0.168 on a 12 3 x 24 lattice. Upper panel: 
ji = 0, "low plaquette"; lower pannel: fi = 0, "high plaquette". 
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Figure 13: Eigenvalues of the Wilson-Dirac fermion matrix with small 
absolute value in case of the Wilson plaquette action at [3 = 5.20, \i = 
0.01, k = 0.1715 on a 12 3 x 24 lattice. Both "high plaquette" and "low 
plaquette" spectra are shown. 
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Figure 14: The schematic view of the phase transitions in the ((3, k, /i) 
space for Wilson quarks with both DBW2 and Wilson plaquette gauge ac- 
tion (f3=bare gauge coupling, K=hopping parameter, fi=bare twisted quark 
mass, fi K = (2/t) -1 =bare untwisted quark mass.) The crosses mark the 
second order boundary line of the first order phase transition surface. At 
strong gauge coupling there is the surface containing the Aoki phase, which 
ends at a point denoted by (5 — The figure does not extend down to 
(3 = and only one "finger" of the Aoki phase is shown. 
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